
#include "bups.h"

#include <cassert>

/*
 * Clase auxiliar que maneja una clase especial de bases conocidas como
 * Bases Uniforme Periodica Spline (BUPS)
 *
 * Los puntos nudo (o knots) de una BSpline son caracterizados por la secuencia
 * de escalares creciente (no necesariamente monotona) t1, t2, ... ti
 *
 * Dicha secuencia en una BUPS cumple que tj - ti = j - i (Uniforme equiespaciada)
 *
 * */

float BUPS::base(int ordenK, int baseI, float t) const {
   if(ordenK == 1) {
      return ((t >= baseI) and (t < baseI+1)) ? 1 : 0;
   }

   float coefI = (t - baseI) / (baseI+ordenK-1 - baseI);
   float coefSigI = (baseI+ordenK - t) / (baseI+ordenK - (baseI+1));

   return coefI*base(ordenK-1, baseI, t) + coefSigI*base(ordenK-1, baseI+1, t);
}

float BUPS::baseDerivada(int ordenK, int baseI, float t) const {
   if(ordenK == 1) {
      return 0;
   }

   float coefI = (t - baseI) / (baseI+ordenK-1 - baseI);
   float coefSigI = (baseI+ordenK - t) / (baseI+ordenK - (baseI+1));

   float coefIDerivado = 1.0 / (baseI+ordenK-1 - baseI);
   float coefSigIDerivado = -1.0 / (baseI+ordenK - (baseI+1));

   return coefIDerivado * base(ordenK-1, baseI, t) + 
      coefI * baseDerivada(ordenK-1, baseI, t) +
      coefSigIDerivado * base(ordenK-1, baseI+1, t) +
      coefSigI * baseDerivada(ordenK-1, baseI+1, t);
}

BUPS::BUPS(int orden) : orden(orden) {
   assert(orden > 0);
}

// El valor t debe estar entre [orden-1; M-orden]
// para algun "M" tal que "M" > "orden"
//
// Cualquier otro valor para t dara resultados inconsistentes.
float BUPS::base(int baseI, float t) const {
   assert(int(t) >= 0);

   float resultado = base(orden, baseI, t);

   assert(resultado >= 0);
   assert(resultado <= 1);

   return resultado;
}

float BUPS::derivada(int baseI, float t) const {
   assert(int(t) >= 0);

   float resultado = baseDerivada(orden, baseI, t);

   return resultado;
}


